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Abstract 

A single neuron is known to generate almost identical spike trains when the same fluctuating 
input is repeatedly applied. Here, we study the reliability of spike firing in a pulse-coupled network 
of oscillator neurons receiving fluctuating inputs. We can study the precise responses of the network 
as synchronization between uncoupled copies of the network by a common noisy input. To study 
the noise-induced synchronization between networks, we derive a self-consistent equation for the 
distribution of spike-time differences between the networks. Solving this equation, we elucidate 
how the spike precision changes as a function of the coupling strength. 
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Reliable information processing requires a code that can be represented and transmitted 
reliably within the precision of the devices. In the brain, single neurons can generate highly 
precise spike trains when they are repeatedly activated by the same fluctuating input l|. 
Neurons, however, work collectively rather than individually in their network. It remains 
unclear whether precise elements put in a network can still respond precisely since mutual 
couplings may affect the response of the individual neurons [2J. To answer this question, 
we investigate the temporal precision of responses of a pulse-coupled network of oscillators 
when a set of independent fluctuating inputs, i.e., frozen noise, is repeatedly applied to the 
network. To study this problem analytically, we introduce uncoupled copies of the network 
that commonly receive the set of fluctuating inputs. Suppose that the original network 
repeatedly generates identical precisely-timed spike responses. Then, we can interpret the 
collection of these responses across trials as responses of the individual copies in a trial. Thus, 
in-phase synchronization between identical networks implies precisely-timed responses of a 
single network across trials. 

Such noise-induced synchronization was previously studied between single oscillators 
0, Here, we study the noise- induced synchronization between networks of oscillators 
to clarify whether each network is able to encode information about fluctuating inputs 
into precisely-timed spikes. Noise-induced synchronization can appear in lasers^], chemical 
reactions J?}], gene networks^, electronic circuits and neural systems. For instance, noise- 
induced synchronization of neural oscillators may play an active role in olfactory information 
processing 9j . An attempt was made to employ such an oscillating synchronization for a 
dynamic clock in digital VLSI devices jl 01] . Thus, clarifying the underlying mechanism has 
significant implications for a variety of dynamical systems. 

We develop a mean-field theory of noise-induced synchronization between the copies of 
the pulse-coupled oscillator network. We can analytically derive a self-consistent equation 
for the distribution of phase differences between the corresponding oscillators and obtain 
the distribution as a function of the coupling strength. In so doing, we assume that the 
average of the connections projecting to each oscillator vanishes. The distribution allows 
us to reveal the nontrivial effect of mutual couplings on the noise-induced synchronization 
of oscillator networks or, equivalently, the temporal precision of responses in an oscillator 
network. 

Let us consider multiple trials in which a pulse-coupled network of homogeneous iV os- 
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cillating neurons receive the same fluctuating input in all trials: 
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where i = 1, • • • , iV, xf 1 ^ is the state variable of the ith neuron in the ath trial, F ^a^jis 
the intrinsic dynamics of the neuron, is the coupling strength from the jth to the ith 
neuron, is the nth spike time of the jth neuron in the ath trial, and (t) is the fluctuating 
input to the ith neuron. F (x) has a stable limit-cycle solution Xq (t) satisfying Xq = F (xq) 
with period T. Isolated neurons thus fire regularly with a firing rate r = 1/T. We use 
the zero mean white Gaussian noise with variance D as the fluctuating input £j (t) , which 
is independent among neurons whereas it should be the same across trials, = 0, 

(£i(t)£j{s)) = D5ij5(t — s). Since Eq. (JTJ) is a stochastic differential equation, we have to 
clarify its interpretation. We use the Stratonovich interpretation, namely, we define the 



white noise as the limit of colored noise with infinitesimal correlation time[ll|]. Components 
of the zero mean matrix gij take both positive and negative values independently, while they 
should also be exactly the same across trials. We denote the variance of the matrix as G, 
{gij) = 0, (gijgki) = G8ik$ji- Even though the same input and the same network are shared 
by all trials, it is unclear whether spike times would coincide across trials because initial 
values are different among trials. 

In order to obtain a unified description of the problem, we apply the standard phase 
reduction methodfl^] to Eq. ([!]). Regarding both the fluctuating inputs and the interac- 
tions as perturbation pf^ to the limit-cycle oscillators, we obtain the following stochastic 
differential equations for phase variables: 

<T = 1 + Z (4 a) ) (t) , (2) 

N 
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where the phase (pf 6 is defined to increase by T for every cycle of x 1 ^ around the limit 
cycle. Natural angular velocity is thus equal to 1 and the spike time satisfies the relation 
<f>i^ (t^J = nT. The phase response function, or the phase sensitivity, Z(<f>) = V a ,0| 3;=£CO ( 9 i), 
quantifies the phase response to perturbations 12|] . 



Spike time difference across trials is quantified as the distribution function of phase dif- 
ferences across trials in the phase description. If corresponding oscillators in different trials 



synchronize with each other in phase, the distribution will be the delta-function and spike 
trains will be the same across trials. If, in contrast, oscillators across trials do not syn- 
chronize perfectly, the distribution will have a finite width, which characterizes variation of 
spike times across trials. The synchronization across trials is an extension of common-noise 
induced synchronization between uncoupled oscillators. It is noteworthy that the synchro- 
nization mentioned here is different from widely studied synchronization within a network; 
rather this is synchronization across trials that are mutually uncoupled by definition. Actu- 
ally, oscillators within the network tend to desynchronize each other because of independent 
fluctuating inputs they received. 

To obtain the distribution function of phase differences, we first derive the Fokker-Planck 
e qU at 10 n S at lsfl ed by the distrihutioo fi Q. Without ,o SS of geoerahty, we foc US on 
phase differences between two trials, ipi = tpf^ — <f)f~\ The Fokker-Planck equation for the 
distribution P (ipi,t) is expressed as ^ = — -^p {A\P) + (A 2 P), where the 1st and the 
2nd Kramers-Moyal coefficient, A\ and A 2 , are defined as moments of ipi increment in an 
infinitesimal time step, A n (ipi) = lim^o ^ ((^i (t) ~ (0))"") L.(n)=7/>,-- We now invoke the 
averaging assumption to calculate the coefficients 13] 5, ll2j, ll4j. With sufficiently small 
perturbations, the distribution function P (ipi, t) varies slowly compared with the oscillator 
natural period, T(= 1/r). We can thus replace the infinitesimal time-step increment in A n 
by the average increment in the period: 

An {A) * ^ ((^ CO - h (0)) n ) |* (0)= *. (3) 

The increment is calculated from Eq. (T5]) as follows. Integration of Eq. (j2J) from the initial 
phase <pf (0) gives 

^ (t) = 4 a) (0) + t + J * dxZ (*)) pf (x) . (4) 

To evaluate (pf 1 ' (x) which appears in the integrand, we substitute Eq. (T4]) recursively into 
the right hand side and obtain 

^ (t) = ^ (0) + t + f dxZ (0) + x) pf (*) (5) 
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up to the 1st order of p\ a ^ 13]. Using the assumption (g^) = 0, we find 

< Pi (x) pf ] (y) > 

= (& {x) o (y)> + ( E ^•^ <J ( x - '£) 5 ( y _ t 
= ^ (* - v) + ( E 4 5 (* - *S) 5 - ) 

\j,n,m I 

= D5(x-y) + GNrS (x - y) , 

where the normalized cross-correlation function between spike trains of trials is defined as 
S (x — y) = \ nm^( x ~~ tln)$(y ~ 4fm)/(-^ r ) / > which yields the distribution of spike time 



difference between trials. Substituting both ipi 



<(2) 



and Eq. (jSJ) into Eq. 



and 



using the above equation, we obtaine coefficients explicitly as A\ (ip) = and 



A 2 ("0) = — — / dxZ (x^ 2 
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From Eq. ([6]), we obtain the explicit form of the Fokker-Planck equation. The distribution 
P is given as the stationary solution of the equation, 



P(ij) = c/A 2 (V) 



(7) 



where c is a normalization constant. 

Because A 2 (ip) still includes the unknown function S, Eq. (I7j) is not a closed form of the 
distribution P. However, we can derive a quite simple relationship between P and S. From 



the relationship 



(,(<*) f t (a) 



nT, we obtain 
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0. Expanding this 



with respect to spike time difference when the difference is small and using the fact that 



the natural velocity of phase is unity, we obtain tf^ — t\ 



,(2) 
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a 1 ) 



ipi, which means 



that the spike time difference is equal to the phase difference. Therefore their respective 
distributions should also be the same; 



P = S. 



(8) 



Using Eq. (jSj) is a type of the mean- 
physics and theoretical neuroscience 
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. ield approximation which has been widely used in 
Unlike the previous studies, however, we do not 
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FIG. 1: (a) (solid lines) The right hand side of the self-consistent equation when q = 0, 1, and 
10, from top to bottom. Dotted line is the left hand side of the equation. Filled and empty circles 
are stable and unstable solutions respectively, (b) Unstable (dashed line) and stable (thick line) 
solutions of the self-consistent equation. 

treat statistics only within the network; rather we can treat statistics between uncoupled 
trials. 

Equations (6-8) give a self-consistent equation for the distribution function P, and for 
S simultaneously. Since P is a T-periodic function, we can use the Fourier expansion to 
simplify the equation. The nth Fourier component f n of a T-periodic function / is defined 
from / (x) = J2n=-oo f n e i2nnx/T . Noting that P is an even function due to symmetry between 
trials, we can obtain the expression P from Eq. (jBJ) and (J7J) as 

P ^ = E^il^| 2 (l-^cos(^))' 
where c is a normalization constant and W n = (1 + qP n ) /(l + qr). The variable q is defined 
as the ratio of variance between internal input and external input per neuron q = NG/D, 
which acts as the control parameter of our problem. We then obtain the self-consistent 
equation for the Fourier component of P as the final form: 

,2, 
JO 

T 

E^Lll^nl^l-WnCOSnV) 

The infinite Fourier series which appears in Eq. (Q and ({TO]) is terminated in finite since 
\Z n \ vanishes generally when n is large. Therefore, Eq. ([TO]) gives a well-defined equation 
for a finite set of variables corresponding to \Z n \ > 0. Putting the solution of the equation 
to Eq. 0, we obtain the final expression of P, and hence S. 

We now examine the above result using the simplest example where T = 2ir and Z {(f)) = 
1 — cos 0. This phase response function is obtained directly from the quadratic integrate- and- 
fire neuron, which is also known as the theta-neuron or the Ermentrout-Kopell canonical 
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FIG. 2: (a) Width of P (^) as a function of the coupling strength G. The width is calculated from 
both theoretically (solid line) and numerically (filled circles), (b) P (ifi) calculated numerically 
when G = 0.0, 0.2, 0.5, and 1.0 (xlO -3 ), from left to right. Solid lines (online red) are theoretical 
predictions, Eq. ([9]). 



for n > 1 to Eq. ( jlOl) . we obtain the self- 

1 - 



model 16j. Putting \Z\\ = 1/2 and \Z n \ 

consistent equation for a variable Pi as Pi = ^1 — a/1 — W^j / (2ttWi) = Hi (Pi). To ensure 
P > 0, the condition Pi < r is required. We plot each sides of the self-consistent equation 
in Figure 1 for various values of q. When q = 0, there is no coupling in the network, the two 
lines have one intersection at Pi = r. Inserting the solution to Eq. ([9]) gives P (ip) = S (ijj), 
which exactly recover the previous result that single neurons always response precisely and 
that single oscillators always synchronize in-phase by common noise|3(]. Note that the Fourier 
series of P is infinite while only the 1st component appears in the self-consistent equation. 
When q > 0, i.e. there are finite couplings in the network, the solution Pi = r turns out 
to be unstable and another stable solution appears in P x < r. Inserting the solution to 
Eq. ([9]), we obtain the distribution with a finite width. As we increase q, i.e. increase 
coupling strength G, the width increases monotonically. We plot the evolution of the width 
in Figure 2a. In parallel to the analysis, we simulate Eq. ([2]) directly and plot the result 
also in Figure 2a. Numerical results are well fitted by the theoretical curve. Figure 2b 
shows distribution functions P obtained from both numerical and theoretical calculations. 
Theoretical predictions agree fairly well with numerical results. 

Positive width obtained here implies that networks of oscillators, in contrast to a single 
oscillator, do not synchronize perfectly across trials even though every trial is driven by 
the same input. Even if all the constituents, i.e. single oscillators here, are faithful to the 
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FIG. 3: Distribution P (ip) of Hodgkin-Huxley neuron with mean current strength / = 
6.5 mA/cm 2 . calculated numerically from Eq. ([1]) for G = (a), 10~ 2 (b), and 10 _1 (c). 

Solid lines (online red) are theoretical predictions. 

input, it is not the case for the system as a whole. Spike sequence generated by the network 
is thus not perfectly the same across trials. Instead, the result implies enough coherence 
of spike trains between different trials, as indicated by the huge peak of the distribution 
around ip = 0. The width of the peak measures a degree of the coherence, up to which 
spike trains can be used reliably. Furthermore, the result tells us how the coherence changes 
qualitatively as a function of the coupling strength in the network. 

If the network encodes information with precise spikes with temporal precision Ai/j and 
that the information is decoded by a downstream neuron innervated by K neurons in the net- 
work, reliability for the decoder neuron is roughly estimated as (P (0) Aip) K , which rapidly 
goes to zero as K increases, except if G is close to zero. Therefore the result infers that 
population codes based on precise spike times require weak or sparse mutual couplings in 
the network. 



As the second example, we consider a coupled network of the Hodgkin-Huxley neurons 171 ]. 
We numerically integrate the system, Eq. (JTJ, and compare the results predicted by theory. 
Since the Fourier component Z n of the neuron decays rapidly to zero when n > 5, we use 
only the first 5 components to solve Eq. ffTUl) . To realize the Stratonovich situation in the 
numerical simulation, we use the Ornstein-Uhlenbech process with a correlation time r = 1, 
which is sufficiently small comparing to the period of the oscillation. Figure 3 shows the 
result. As similar to the previous example, the distribution starts from the delta-function 
when G = and grows to a distribution with a finite width. The theory agrees well with 
numerical results, especially when the width of the distribution is small. 

We have assumed the balanced couplings, i.e., (g^) = 0. If the couplings are not balanced, 
coefficient A 2 acquire additional terms, which in turn modify the self-consistent equation. 
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Such biased couplings may significantly change the intra-network synchronization property 



of oscillators 



18| . which presumably influences the spike precision or the synchronization 



across trials. It is intriguing to extend the present analysis to the case with (g^) ^ 0. 

Since the self-consistent equation has a unique stable solution, the obtained distribution 
is always realized over repeated trials even if initial phase differences between trials are 
infinitesimally small. This property reminds us of chaos theory in dynamical systems, in 
which a small difference between orbits increases rapidly until it roughly converges to a 
finite value characterized by the size of chaotic attractors. Recent studies sug gested that 
neuronal circuits can be chaotic in a manner useful for computations 19j, |2C|, l2l|- It remains 
fascinating whether our results reveal an active role of chaos in the brain. 
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